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An H theorem for contact forces in granular materials 
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A maximum entropy theorem is developed and tested for granular contact forces. Althougli it is 
idealized, describing two dimensional packings of round, rigid, frictionless, cohesionless disks with 
coordination number Z = 4, it appears to describe a central part of the physics present in the more 
general cases. The theorem does not make the strong claims of Edwards' hypothesis, nor does it rely 
upon Edwards' hypothesis at any point. Instead, it begins solely from the physical assumption that 
closed loops of grains are unable to impose strong force correlations around the loop. This statement 
is shown to be a generalization of Boltzmann's Assumption of Molecular Chaos (his stosszahlansatz) , 
allowing for the extra symmetries of granular stress propagation compared to the more limited 
symmetries of momentum propagation in a thermodynamic system. The theorem that follows from 
this is similar to Boltzmann's H theorem and is presented as an alternative to Edwards' hypothesis 
for explaining some granular phenomena. It identifies a very interesting feature of granular packings: 
if the generalized stosszahlansatz is correct, then the bulk of homogeneous granular packings must 
satisfy a maximum entropy condition simply by virtue of being stable, without any exploration of 
phase space required. This leads to an independent derivation of the contact force statistics, and 
these predictions have been compared to numerical simulation data in the isotropic case. The good 
agreement implies that the generalized stosszahlansatz is indeed accurate at least for the isotropic 
state of the idealized case studied here, and that it is the reductionist explanation for contact force 
statistics in this case. 

PACS numbers: 45.70.Cc, 05.20.Gg, 05.65.+b 



I. INTRODUCTION 



Edwards hypothesized that every mechanically stable 
micro state of a powder is equally probable in a Gibbs- 
like ensemble, for cases where the powder is prepared in 
certain, repeatable ways 0]. This hypothesis has been 
extended by others to the statistics of granular contact 
forces for both isostatic and hyperstatic cases. Packings 
of round, rigid, frictionless disks or spheres are known 
to be isostatic [2j], meaning that there are exactly the 
same number of contact force unknowns in a granular 
packing as there are stability equations for the grains, so 
that the value of each contact force can be resolved by 
linear algebra operating solely upon the geometry of the 
inter-granular contact network (assuming that the over- 
all stress of the packing has also been specified). For 
isostatic cases, then, each valid micro state in an Ed- 
wards' ensemble corresponds to exactly one micro state 
of contact forces. An ensemble of packings with a flat 
measure thus implies an ensemble of contact force net- 
works with the same flat measure. Hyperstatic packings, 
on the other hand, have more contact forces unknowns 
in the granular packing than there are stability equa- 
tions for the individual grains. Thus, the values of the 
forces throughout the packing are mechanically indeter- 
minate. A Gibbs-like ensemble method has been applied 
to the contact forces in such packings by selecting a single 
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packing geometry that is held constant throughout the 
ensemble, assuming a flat measure over all of the valid 
contact force microstates [1, All these ensembles, iso- 
static as well as hyperstatic, are similar to Gibbs' ensem- 
bles in thermal systems because the flat measure is as- 
sumed a priori. The isostatic case was recently analyzed 
by extending Gibbs' most probable distribution method 
d, 0] , the same method used in textbooks to derive the 
Maxwell-Boltzmann, Bose and Fermi distributions. In 
those derivations, the most probable distribution is con- 
cerned with momenta /(p) or particle energies f{E), but 
when extended to granular contact forces it predicts the 
most probable distribution of single grain states, p{g), 
where 5 is a set of variables that describe everything 
that can be known about an individual grain including 
all of its contact angles, forces and their correlations. 
Subsequent discrete element modeling has validated the 
predictions of this theory for the special, isostatic case 
described above [1,01 • Thus, Edwards hypothesis is suf- 
ficient to derive granular contact force statistics for single 
grain states in this case. 

However, this paper will show that Edwards' hypoth- 
esis is not really necessary and that it is not really a 
part of the reductionist explanation for granular contact 
force statistics. It makes a very strong claim about mi- 
crostructural geometries being equally probable, and this 
is a much stronger claim than we need to derive granular 
contact force statistics. Instead, the reductionist expla- 
nation is found in the nature of correlations in granular 
contact forces, which can be summarized in a statement 
that is a generalization of Boltzmann's Assumption of 
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FIG. 1: Illustration how sheets of force are conserved while 
translating through a 2D granular packing (as described in 
the text), so that a Boltzmann-like transport equation may 
be defined. 

Molecular Chaos. Boltzmann showed that the absence 
of pre-correlation between colliding gas molecules inex- 
orably relaxes the gas to maximum entropy; this paper 
will show that the corresponding condition in granular 
contact forces inexorably requires the packing to exist 
(in its bulk) in a relaxed state of maximum contact force 
entropy. This paper will follow Boltzmann's approach, 
developing a "translation" equation analogous to Boltz- 
mann's transport equation (but adapted for static pack- 
ings in which nothing is moving), with a maximum en- 
tropy proof similar to the H theorem, to produce a new 
derivation of p(g) without making any a priori assump- 
tions of a flat measure. This new theory is presented as 
an alternative to Edwards' hypothesis to provide a more 
physical basis for the statistical mechanics of some gran- 
ular phenomena. 



II. THE GRANULAR CONTACT FORCE 
S TOSSZA HLA NSA TZ 

A. "Translation" in a Static Granular Packing 

Boltzmann wrote a transport equation to track the 
evolution through time in the single particle distribution 
function f{p), where p is momentum. For static granu- 
lar packings we must describe how the statistics evolve 
through space, not time. We will therefore perform 
"translation" through successive cross-sectional "layers" 
of grains, as defined here. Referring to Fig. ([T]), we: (1) 
draw a cross-sectional line (like any one of the straight 
lines in the figure), (2) select the set of grains intersected 
by that one line (shaded grains), and (3) include only 
those contacts on that set of grains (heavy dots) that 
(a) connect to other grains external to the shaded set 
and (b) are located on one side of the shaded set. Next, 
we consider the force vectors on that set of contacts, in- 
cluding both the normal and tangential components of 
force if the grains are frictional. We decompose these 
vectors into Cartesian components parallel and normal 
to the dashed line. Summing all the components in the 



normal direction produces the "total Cartesian load" in 
that plane. "Translating" refers to the continuous mo- 
tion of the cross-sectional line that picks out the set of 
grains in a layer. As the line translates across the pack- 
ing, some grains are no longer intersected by the line and 
hence leave the layer, while some other grains become 
intersected by the line and hence join the layer. In the 
absence of gravity, the total Cartesian load is conserved 
with respect to translating the straight line, as illustrated 
in the figure by the parallel layer of shaded grains and 
heavy dots. If that were not so, then the grains between 
the two shaded layers would be accelerating. 

Because the container walls might also contribute 
forces parallel to the direction of translation, thus spoil- 
ing the conservation of the total force, we must con- 
sider special cases where this cannot occur. There are 
at least two such cases, corresponding to microcanoni- 
cal and canonical statistics. If the container is 2D (for 
specificity), and if the side walls are flat, parallel and fric- 
tionless, then we will have two spatial axes along which 
to translate layers of grains subject to the conservation 
laws as shown in the left and right sides of Fig. [TJ We 
cannot perform the translation in a diagonal direction 
in this finite container because then the length of layers 
will not be constant and neither will the total perpendic- 
ular force contained within them. This translation in a 
finite container corresponds to microcanonical statistics 
because the system is closed and the force conservation is 
exact. Alternatively, we can perform the translation in an 
extremely large packing (perhaps infinite) and consider 
only the set of grains selected by a relatively short line 
segment translating through the middle of the packing far 
from any boundaries. In that case, the conservation law 
will not be exact because some force will be entering and 
exiting the segment at both of its ends throughout the 
translation. However, the longer the segment of grains 
is, then the less significant those fluctuations become and 
the more closely it approximates an exact conservation 
law. This corresponds to canonical statistics because the 
force contained in the finite line segment is in contact 
with an infinite reservoir of force at each of its ends. In 
this canonical case the line segment may have any ar- 
bitrary orientation, not only the x or y directions, and 
it translates across the packing in a direction normal to 
its orientation. In this paper it does not matter whether 
we consider the microcanonical or canonical case, nor 
in what direction we imagine the translation to occur. 
The generality of the equations describing this transla- 
tion contributes to the generality of the conclusions. 

It should be noted that this translation ansatz allows 
us to treat a static granular packing with methods bor- 
rowed from kinetic theory, even though the packing is 
completely static. The role of the time dimension from 
kinetic theory is replaced by the translation through a 
spatial dimension. The collisions of particles from ki- 
netic theory are replaced by the meeting together of con- 
tact forces on the grains. The force vectors on one hemi- 
sphere of a grain are analogous to gas particles going into 
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TABLE I: Compact notation for 2D packings with Z = A. 
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a collision, and the forces on the opposite hemisphere are 
analogous to the same particles exiting from the collision 
later in time. (The number of contact forces "entering" 
and "exiting" the opposite hemispheres of a grain in this 
fashion need not be conserved.) Thus, the granular the- 
ory developed here is analagous to molecular kinetic the- 
ory, but it must be remembered that everything is com- 
pletely static and that all contacts between grains are 
unchanging in time, not transitory like the collisions in 
kinetic theory. This theory is not to be confused with the 
kinetic theory of granular gases, which deals with grains 
moving and colliding in time. 

Although this force conservation applies to more gen- 
eral cases, the remainder of this paper deals only with 
round, rigid, frictionless, 2D grains. Table [J defines the 
compact notation that will be used to describe the states 
of (1) single grains, (2) sets of grains, and (3) entire pack- 
ings in this special case. The single-grain state variable 
is 

g^{wx,Wy,0i,...,0i) (1) 

where the Cartesian loads Wx and Wy are the total force 
borne by a grain in each orthogonal direction and Oi are 
the contact angles. From these variables, the individual 
contact forces may be recovered by linear algebra and 
trigonometry. This form of g assumes contact number 
Z = 4 for every grain. It may be generalized for grains 
with Z 7^ 4 but doing so adds no insight to the physics 
at this stage. Note also that certain regions in g describe 
grain configurations that have one or more negative (ten- 
sile) forces. It will therefore be necessary to restrict the 
range of g to the non-tensile (stable) region S when deal- 
ing with cohesionless grains, as we do in this paper. 

The densities of states in Table U refer to either single 
particle or multiple particle states in the corresponding 
phase space. For example, the density of single particle 
states p{g) refers to the grains in a single packing or in 
a single layer of a packing (depending on the context), 
and it tells how many of those grains exist per unit vol- 
ume of g-space as a function of the location in g-space. 
Likewise, the density of packing states p{T) refers to a 
statistical ensemble of granular packings and tells how 
many of those packings exist per unit volume of F-space 
as a function of the location in F-space. The construc- 
tion list C contained in F specifies the exact ordering in 
which the grains are connected together such that they 
form a packing. The density of states /o(7) is also a multi- 
particle density of states describing collections of grains 
with specified single-particle states, but unlike p(F) it 




Integrate across _ _ _ _ 

time and point all /, + /, + + /a = 

momentum vectors 

inward Identical to forces on 

static grains 

FIG. 2: Illustration how momentum vectors in binary col- 
lisions of a dilute gas are analogous to contact forces on a 
static grain. Boltzmann wrote the stosszahlansatz with pref- 
erence for the in-going momentum vectors, but for static gran- 
ular materials the symmetries among the contact force vectors 
must be restored. 

does not have a construction list and so it does not specify 
the location of the grains in physical space nor whether 
they physically connect to one another in a packing. De- 
spite the mathematical abstraction, it will be useful in 
this paper. 

B. The Need to Generalize Boltzmann's 

Stosszahlansatz 

Boltzmann's stosszahlansatz, or his Assumption of 
Molecular Chaos, is that colliding particles have uncor- 
related momenta prior to the collision, 

FiPl,P2) = fiPi) ■ f{P2) (2) 
or more compactly, 

F21 = /1/2 (3) 

where F21 is the joint probability distribution for the 
two particles and fi are the single particle distribution 
functions. Are the contact forces that meet together on a 
grain in a granular packing similarly uncorrelated? Fig. [2] 
shows how the relationship of contact forces upon a static 
grain does bear an analogy to momenta in a gas. In a 
particle collision the sum of the momentum vectors is 
conserved. Hence, if the outgoing momentum vectors are 
reversed, then the sum of all four momentum vectors will 
equal precisely zero. Putting a little space between the 
four arrow heads and drawing a circle, we see that this 
set of four vectors is identical to four contact forces on a 
static grain. 

However, Eq. [5] assumes non-correlation for the in- 
coming momentum vectors only Q, whereas the ideal- 
ized granular packings considered here (without gravity) 



should be more symmetric than this. If we try to apply 
Boltzmann's stosszahlansatz to granular contact forces, 
we could group the four forces on the grain into six differ- 
ent pairs, and there is no reason why one of those pairs 
should be written as uncorrelated to the exclusion of the 
other five. Furthermore, Silbert, Grest and Landry have 
demonstrated through numerical modeling that contact 
forces meeting together on a grain do have a very strong 
pattern of correlation and anticorrelation There is 
anti-correlation for contacts closer together than roughly 
7r/2 radians of angular separation, and positive correla- 
tion when the angular separation is greater than roughly 
7r/2. 

Because of these things, it would be incorrect to use 
the form of Eq. [5] for the granular stosszahlansatz. In- 
stead, we must begin with a statement that is even more 
fundamental than Eq. [2l namely, that correlation can- 
not arise through the closure of loops of forces (or loops 
of momenta) within the network of these vectors. This 
is the generalized stosszahlansatz. It will be shown be- 
low that this generalized statement does indeed reduce 
to Eq. [2] when causality is assumed to proceed only in 
the forward time direction. However, it produces a very 
different form when "causality" (information flow) is as- 
sumed to be symmetric in all dimensions of the vector 
network. That symmetric form is the one that we need 
for granular packings. 

A priori arguments have been provided in [5] to jus- 
tify the generalized stosszahlansatz, that correlation does 
not arise through the closure of force loops in the force 
vector network. Appendix A of this paper extends those 
arguments. As with thermal systems, the ultimate proof 
of the stosszahlansatz will be its ability to make correct 
predictions. If the predictions are correct, then we may 
claim a posteriori that it is indeed the reductionist ex- 
planation for the statistics of granular contact forces. 



C. Mathematical Form of the Granular 

Stosszahlansatz 

For a dilute gas it is possible to represent the net- 
work of particle collisions occurring through time as a 
4-dimensional network graph, with the collisions as the 
vertices and the particle trajectories as the edges (line 
segments) connecting the vertices, and with one of the 
four graphical dimensions representing time. A lower- 
dimensional version of this is illustrated in Fig. O It will 
not be shown here, but the density of states for a block of 
time in this network can be written using variables that 
are identical to those used for granular packing densities 
of states, with the distribution of momenta f{p) replac- 
ing the distribution of contact forces P{f ) and the time 
dimension replacing any one spatial dimension, and with 
some differences in the term that defines the allowable 
sets of momentum vectors (or force vectors) that appear 
at the collisions (or at the grains). It is possible to an- 
alyze the ensemble of these networks in each case and 
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FIG. 3; Illustration how a network of gas collisions and a 
network of packed grains may both be represented by network 
graphs of similar topology. Identical sets of variables may be 
used to define a density of states of the vector networks for 
either system. 



compare the results. Here we will show how two different 
forms of the stosszahlansatz are obtained depending only 
upon the directionality of the information flow through 
the network. We will obtain Boltzmann's stosszahlansatz 
when we assume all information travels in one direction 
across the network (representing time-asymmetry), and 
we will obtain the granular version of the stosszahlansatz 
when we assume information has propagated across the 
network symmetrically in all directions. This will demon- 
strate that the generalized stosszahlansatz really does en- 
compass both versions and reduces to one or the other 
depending on the assumed direction of information flow. 
The density of states we wrote in this section will then 
no longer be needed in the remainder of the paper. It 
will only be used in this section to derive the two math- 
ematical forms of the stosszahlansatz. 

The multi-particle density of states describing the en- 
semble of these granular packings (and adaptable to the 
collision networks of dilute gases) will be written in a 
form like this. 



Pi^\{Pj}) = (^(Newton's Third Law) 
x9(No Tensile Forces) 
X 8 (Geometric Closure of Loops) 
X(5(Sequence {Pj}) 
X (5 (Collision Auxiliary) (4) 



For static granular packings of rigid, round, frictionless 
grains (without the symmetry- breaking effect of gravity). 
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this becomes 

PiT\{p,{f)},p,o) = I n Kr^ + f'')] 

X |nn0(r^)| 0c(ric) 

ta=l£=l J 

X X{5{m)--PAm) 

X 5{P^e - Qi0(T)) 

(5) 

Most of the terms in this equation are just notional and 
do not have to be defined in enough detail for doing 
calculations because they will factor out of the equa- 
tions below and will not affect the results in this pa- 
per. Before describing these terms, we note that Eq. [5] 
does not have a term to enforce Newton's Second Law 
on each grain. Instead, we have defined g and T such 
that Newton's Second Law is automatically satisfied for 
every grain throughout the range of the phase space. To 
do this, g was defined with Wx and Wy rather than four 
contact forces. This removes two degrees of freedom per 
grain from the phase space coordinates and thus pre- 
vents those combinations of contact forces on a grain 
that would not have summed to zero. The individual 
contact forces may be calculated for any grain by linear 
algebra and trigonometry from the Wx and Wy with the 
four contact angles contained in g Thus, the con- 

tact forces appear in several terms of Eq. [5] as functions 
^at = f°"^{ga) although the argument was suppressed 
for compactness. 

The first term in Eq. [5] ensures that the density of 
states is nonzero only in the regions of F where Newton's 
Third Law is satisfied for every contacting pair of grains. 
The contact force on grain a on its contact e is f°"^. 
The construction list C tells which grain/contact (ae) 
connects to grain/contact (/3(S) and thereby constructs a 
specific packing from the set of grains specified in F. This 
is the most important term affecting the results of this 
paper. 

For the second term, the Heaviside step function en- 
sures that the density of states is nonzero only where the 
forces are positive, with positive defined as pointing in- 
wardly on the grains. This is because we are dealing with 
cohesionless granular materials (no tensile forces possi- 
ble). 

For the third term, the notional function 0^(F|C) was 
introduced by Edwards 1]. It was defined to evaluate 
either to zero or unity, being unity only when the locus 
in F describes a configuration of grains that form pre- 
cisely closed loops, with contacting grains just touching 
each other precisely without any overlapping. This, along 
with the constraints on the contact forces, is sufficient to 
ensure a physically stable granular packing. has not 
been defined more specifically because it is not needed in 
this paper. 



Discussion of the fourth term will be delayed until last 
due to its importance. The fifth term is the one that 
was called the "collision auxiliary" in Eq. [D This term 
will take different forms depending on whether we are 
discussing the case of granular packings (contact force 
networks) or the case of dilute gases (momentum net- 
works). The similarities and differences between the two 
cases are unimportant to the theory developed in this pa- 
per, but they are explained in the Appendix if the reader 
is interested in the analogy between the two cases. It 
was called the "collision auxiliary" because in the case 
of dilute gases it would provide more constraints on the 
set of physically realizable collisions. In Eq. \5\ the term 
is written in a form relevant only to granular packings. 
The delta function ensures that the density of states is 
nonzero only for packings that have the fabric distribu- 
tion P/^g specified within its argument. This fabric is the 
joint contact angle distribution P4e(^i, ^2, ^3, ^4) provid- 
ing the probability that all four contacts on an individ- 
ual grain will take on specified values simultaneously, as 
discussed in This joint distribution is important 

because of the intra-grain correlations between contact 
forces. The function Q40 computes this statistical distri- 
bution from the set of grains {ga} specified by the locus 
F. In the thermodynamic limit with an infinite number of 
grains, the delta function selects only those packings that 
have precisely the specified P^e- For a practical packing 
with only a finite number of grains, this term should be 
defined to allow some statistical fluctuation in P^g . How- 
ever, as stated above, we do not need to define this term 
so precisely since it is just a notional placeholder and 
does not affect the rest of the paper. 

The fourth term may be the most confusing, and it 
factors out of the equations and does not affect the re- 
sults of this paper, and yet it is the most important term 
for understanding the purpose of the ensemble. Why do 
we begin a paper on a Boltzmann-type transport equa- 
tion with an ensemble, which is reminiscent of Gibbsian 
methods? We do so because we will apply the gener- 
alized stosszahlansatz to this ensemble and then ana- 
lyze its statistics to obtain a mathematical form for that 
stosszahlansatz. That mathematical form can then be 
used in the Boltzmann-type transport equation and we 
can abandon the ensemble methods at that time. But for 
now, this fourth term defines a set of packings in which 
the density of single particle states can evolve with re- 
spect to translation in one of the spatial dimensions, just 
as Boltzmann was concerned with dilute gases in which 
the distribution of momenta would evolve with respect to 
translation through the time dimension. Therefore, we 
subdivide the packing into regions that are subscripted 
by i so that we can translate through these regions by 
incrementing j. These regions are defined between suc- 
cessive cross-sectional planes cutting across the packing 
as shown in Fig.[l] The delta function ensures that 
the density of states is nonzero only for packings that 
have the specified distribution of contact forces Pj{f) 
within each region. At this stage of the paper we do not 
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care what the forms of these specified distributions Pj ac- 
tually are or how they may evolve. All we care about is 
that every packing in the ensemble has the same sequence 
of distributions {Pi{f), Piif), P3{f)i • ■ •} as all the other 
packings in the ensemble. Later this paper will show how 
the sequence of Pj actually does behave 13]. The func- 
tions 'Pj{f\r) compute the distribution of contact forces 
from the grains {ga} in each region j specified by the lo- 
cus r. As with the similar function Q^g discussed above, 
for practical packings with only a finite number of grains 
the term should be defined to allow some statistical fluc- 
tuation in the contact force distributions. However, this 
term is notional and will factor out from the equations 
below. 

Finally, we note that there is no term to specify the 
stress state of the packing. That is because it is specified 
implicitly in the {Pj}- We can compute the stress in any 
region as a function of the contact force vectors in that 
region. Newton's third law enforced across the packing 
will ensure that any stresses appearing in adjacent re- 
gions are physically possible. 

Because the identities of the grains are unimportant, 
the statistics of the ensemble are unaffected by summing 
the density of states p over any set of permutations of 
grain exchanges. At any particular location in F-space, 
the delta function that enforces Newton's third law will 
select from those permutations only those that produce 
a self-consistent packing. Here for convenience we sum 
over those permutations that keep the grains within their 
own regions j, denoted by the subset {i}* of the set {i} 
of all possible construction lists Ci. 

p{T) ^ J2 /'(r) (6) 

All terms in Eq. [5] factor out from this sum except those 
containing Ci explicitly. We shall examine this unfactored 
part in the thermodynamic limit 

'^' = ,.,1™ ^ E 0c(r|c.)n'^'(/"" + /'') (7) 

{N,}-,{^} ^^^^^^ C. 

Now we apply the generalized stosszahlansatz by re- 
moving the closure of force loops in the packings. We 
do so simply by eliminating the term from Eq. [5l 
as explained by Fig. H) This converts the ensemble of 
closed-loop contact force networks into an ensemble of 
tree networks as shown in Fig. |4l If the generalized 
stosszahlansatz is correct, then eliminating will have 
no effect on the statistics of the contact forces or the 
single grain states in the ensemble. 

We may note that it was a strategic choice to specify 
the joint distribution of contact angles P^e as the fabric in 
the fifth term, because it implicitly enforces steric exclu- 
sion by evaluating to zero everywhere that adjacent grain 
contact angles are too close together. Therefore, while 
0^ served the purpose of enforcing the non-overlapping 
of grains in the packing, P^g continues to serve this pur- 
pose within the first coordination shell of every grain. 



No closure of force loops 




a b 



FIG. 4: (Left) When 0j is kept as part of the density of states, 
forces form loops and so intersecting forces on the central 
grain may have pre-correlation, because they have seen parts 
of each other before. (Right) When is omitted from the 
density of states, then force loops do not close and so the 
forces that meet together on a grain cannot be pre-correlated 
by the packing. The only correlation arises from the stability 
requirements of the central grain, itself. 




0g enforces steric exclusion in P^q {9i, 62, 63, ©4) enforces 
the second coordination shell steric exclusion in the first 

and higher coordination shell 

FIG. 5: Illustration of coordination shells. The central grain 
is shaded dark. The first coordination shell is shaded lightly. 
The second coordination shell is dashed and unshaded. Omit- 
ting Q( from the density of states removes the closure of force 
loops but does not affect the fabric of the individual grains. 



as illustrated in Fig. [51 Therefore, even though we have 
eliminated 0^, we still have steric exclusion in the first 
coordination shell and so the set of individual grains in 
the ensemble will still be valid. Overlapping of grains 
will appear only in the second and higher coordination 
shells [3. 

Here it is most important to note that Eq. [S] is not an 
Edwards ensemble, which would have included all possi- 
ble sequences {Pj}- We do not specify the relative prob- 
ability of systems having one particular sequence versus 
any other; we only specify one particular sequence to be 
retained in this subset of the Edwards ensemble. The 
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flat measure within this subset of the Edwards ensemble 
is inherent to Boltzmann's kinetic theory when we as- 
sume his stosszahlansatz. That is, in the thermodynamic 
Umit when we drop the closure of loops every tree net- 
work becomes identical apart from a trivial relabeling of 
the branches. The stronger claims made by Gibbs' and 
Edwards' flat measures have been avoided in favor of this 
much weaker claim inherent to Boltzmann's theory. 

The analysis of $ is shown in the appendix. Two dif- 
ferent results can be obtained depending on the assumed 
direction that information has propagated through the 
network. If we assume that information propagated only 
in the direction of increasing j, then we obtain 

F,{hj2)^Pj{h)-Pj{f2) (8) 

which is Boltzmann's form. We may also derive the sym- 
metric case on the same network, with the information 
propagating symmetrically in both directions of j. This 
situation describes, for example, a linearly elastic system 
of grains after the dynamic stress waves have dissipated 
and the solution to the Laplace equation remains. As 
shown in the appendix, we obtain 

FM = [PAh) ■ PAh) ■ Pjih) ■ PAh)Y'^ (9) 

where fi = fi{g) and where the four forces can be located 
with complete generality around the grain (including 3 
contacts on one side of the grain and only 1 contact on 
the opposite side, etc.). This is the mathematical form 
of the granular stosszahlansatz. This form is now sym- 
metric over all four forces, and yet because of the square 
root it still has dimensions of [/]~^. This is so simple that 
perhaps we should have guessed the result before actually 
deriving it. Clearly we did not need to derive the classical 
form given in Eq. [51 since it is well-known from probabil- 
ity theory that the joint distribution of two independent 
variables is the product of their individual distributions. 
However, the symmetric form given in Eq. [Hlis unknown 
in probability theory (to the author's knowledge) because 
probability theory deals with macroscopic events in time, 
in which causality occurs asymmetrically in the forward 
time direction, only. Eq. [5] is the first new result of this 
paper. 

III. GRANULAR "H" THEOREM 
A. The Granular Translation Equation 

We shall analyze what happens to the density of layer 
states defined as in Fig. [T] as its cross sectional line trans- 
lates in some direction x, which is perpendicular to the 
layer and corresponds to increasing j for the successive 
Pj{f)- Therefore, p{g) = pj{g) = p{g,Xj) although the j 
and the x shall be suppressed. 

As we translate from with 
Aa; << -Dparticie (the grain diameter), a small fraction 
of the grains in the layer will no longer be intersected by 



the line and hence will exit the layer. Also, some new 
grains will be intersected by the cross-sectional line and 
hence they will join the layer. If the layer contains M 
grains, then the number of grains expected to leave the 
layer in Ax is to = MAx/Dparticic- If the fabric is con- 
stant across the packing (so that M is constant), then 
this is also the number of grains joining the layer in Ax. 
For sufficiently small Ax the grains exiting the layer will 
be sufficiently far apart to be statistically independent. 
This avoids the need to explicitly account for correla- 
tions in the layer. The probability for a particular set of 
TO grains to leave during Ax is therefore 

Pout{i)^X{p{go.) (10) 

a=l 

The probability for a particular set of m grains to en- 
ter during Ax can be written in terms of the generalized 
stosszahlansatz from Eq. [9l except that we must be care- 
ful because in general P{f ) for contacts that are on one 
hemisphere of the grains will not be the same as for the 
other hemisphere because p{g) is evolving with x and 
P = P[p]. Therefore, we write the stosszahlansatz in the 
form, 

T^(5) = n ^P'ifnig)) n ^PiMa)) (11) 

J)+ 77- 

where 77^ refers to all the contacts on grain a in the 
reverse direction of the translation, and i]'^ refers to its 
contacts in the forward direction of the translation. Like- 
wise, P is the distribution of contact forces in the reverse 
direction and P' the distribution in the forward direc- 
tion. Because the grains entering the layer are statisti- 
cally independent for sufficiently small Ax, we may write 
the probability for a particular set of to grains to enter 
during that interval as 

m 

PM)=^fX{^^{9'o.) (12) 

a=l 

where M is for normalization. 

Because forces are conserved from one layer to the 
next, the physics are analogous to Boltzmann's binary 
collisions, except that they are "m-nary" transitions of 
grains entering and exiting each layer, and that each iter- 
ation Ax contains only one of these TO-nary transitions. 
The probability of having a particular transition 7—^7' 
is 

m 

n(7-7')=AAnp(5o)T±(5:,) (13) 

a=l 

The probability that 7 will go to any set of single grain 
states is 

«(7^all) = / Vg[---Vg'^ ^f U p{9c.)r^{9L) 

(14) 
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where the integrals are carried out only over the stable 
region 5 in g^. Likewise, 

,,(aU -> 7) = / -Dg', ---Vg'^ [] P(5a)T±(3a) 

■^•S a=l 

= / 'Dg[---Vg'„, N 
Js 

\ a a / 

(15) 

To determine the rate of change in p(g) during the 
translation, we write, 



We define H as 



P(5i) = / 'Dg2---Vg„i [n(all ^ 7) - n(7 ^ all)] 



da; 



which is the "translation equation." 



(16) 



B. Counting States for Granular Packings 

To evaluate how this translation equation behaves we 
need a functional similar to Boltzmann's H. As a math- 
ematical proof, the H theorem does not demand that H 
have any physical meaning. However, for the dilute gas 
H becomes the negative of Shannon's entropy when the 
system is in equilibrium, and it is helpful after the math- 
ematical proof is complete to discuss the physical mean- 
ing of this. Following Boltzmann, the granular H shall 
likewise be (the negative of) a generalization of Shan- 
non's entropy. We define H so that it will indicate how 
many packing states correspond to any particular ^(17). 
The "most entropic" p{g) is the one that arises in the 
greatest number of packings. As explained in Ref. [sj, we 
may explicitly count the packing states as a functional 
of p{g) except that here it shall be applied to a single 
layer instead of an entire packing. First we discretize 
p[g) — > Vijkimm where the six arguments of g have been 
broken into small bins of size Aw and and indexed 
as Wxi^ Wyj, Qik, ■ ■ ■ i Oin- Each bin is further divided into 
s smaller bins to enable the typical binomial counting 
without Pauli exclusion, 



'ijklmn 



}-nnnnnn 



/ m n 



(s - 1 + i/i...„)! 



[s-iy. ily^...ny■ 



[T,...„«', 



(17) 



See ^ for the details, s shall drop out of the equations 
in the continuum and thermodynamic limits when Ag — 
(Aw;)^(A6')'* ^ as — * cxD. evaluates to zero when 
any of the g states imply tensile forces, and it evaluates to 
unity otherwise. For simplicity we will restrict all further 
mathematics to the stable region S so we may drop ^' 
from the expressions. 



lim log il 

N ^ 00 
Ag^O 



(18) 



The natural logarithm of il, using Stirling's approxima- 
tion, may be written 

log n = J2^,j,k,l,m,n " 1 + "i-n) log(s - 1 + Vi,„^) 

-{s ~ l)log(s - 1) - t^i...„logt/i...„ 

+ l'i...nl0gTi...„] (19) 

Expanding the first term in a Taylor series around fi^^n — 
0, setting s = NAg, and taking the continuum and ther- 
modynamic limits such that s » i'i,..n, we obtain 



H 



Vg p{g)\og^ 



(20) 



The functional H = H[p(g)] is a measure of the number 
of states in ^(7) that correspond to p{g), and is in fact 
(the negative of) the generalization of Shannon's entropy 
for granular contact forces. In the translation, p{g) must 
be allowed to evolve layer- by-layer, and therefore so must 
P(/). Hence, we distinguish between the hemispheres of 
a grain and use the form of T from Eq. [11] to write. 



H 



(21) 



C. Behavior of p and H in the Translation 

With this metric in hand, the question we wish to ad- 
dress is whether (d/da;)iJ < during the translation de- 
scribed above. Differentiating H, 

T±(s)- 



ax jg ax 



14- log ■ 



p{g) 



+p{g)—logT{g) 
ax 

The last term (call it x) niay be expanded, 
X = ~ J T^g Pig)^^ogT{g) 



(22) 



-Dg p{g)^^ogl[^p{f,{g)) 

-^E^^5 p{g)^\ogP{U9)) 



(23) 



and then evaluated by changing the integration from a 
sum over all grains to a sum over all contacts. 



X = 



/ df P{f) 
Jo 



l_d_ 

2dx 



(24) 
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Substituting Eqs.[Tl[T5l and [16] into {A/dx)H, 

-H ^ M ( Vgi---Vgra I Vg[...Vg'„^ 
' Js Js 



d 

da;' 



1 + log 



.T^(gi) 



pidl) 



(25) 



Since we integrate over S for all ga and g'^ , we may swap 
variables and average the various equivalent expressions 
to obtain, 



d ^ 
da; 



2m ' 



/ Vgi---Vg,n I Vg[...Vg'^ 
l[p(9<^)^^i9L)-l[p(9Ln^i9a) 



X lot 



n.P(5a)T±(g;,) 



(26) 



By inspection of the integrand we see that it is never 
positive for any part of the region of integration. Hence, 



^H<0 
ax 

Furthermore, (d/da;)_ff = if and only if 



(27) 



l[p{g'JT^{ga) = l[p{gc)r^{g'J V g^,g'^ (28) 

a a 

and this is true if and only if {d/dx)p{g) — for all g. 
When that is the case, then P — P' . This proves that the 
bulk of the packing must exist in a relaxed state. More 
will be said about this state, below. 

Furthermore, as shown in the next section, Eq. [28] de- 
fines a sufficient condition to solve p{g) and so by Eq.l^it 
is also sufficient to evaluate H . This produces the small- 
est possible value of i/, since for any other value of H 
(d/dx)iJ 7^ 0. Since this derivation is valid in every pos- 
sible orientation of the layer (because the derivation was 
general and because the stosszahlansatz is symmetric), 
then it must also be valid in the direction of decreasing 
a;, 



d(-a;) 



H <0 



(29) 



This is not true for the dilute gas, because Boltzmann's 
stosszahlansatz is not symmetric [81] • But for the granular 
case with the greater symmetry. 



d d 
—H < and —H > 

da; da; 



which implies 



da; 



H = 



(30) 



(31) 



must be true always, for the bulk of an infinitely large, 
homogeneous granular packing in the special case con- 
sidered here. No exploration of phase space is required 
for the packing to relax. Relaxation is assured through 
spatial relationships, not temporal ones. 

The reason this theory is limited to the bulk of the 
packing, rather than predicting relaxation of stresses 
moving away from the boundary of a container, is that 
there are two things that can cause H to change: evolv- 
ing stresses, and evolving fabric. At the present we do 
not know how to predict the evolution of fabric in the 
boundary layer, and so it is impossible to separate out 
the effect of the relaxation of stresses in that same re- 
gion. By assuming that the fabric is constant as in the 
bulk of a very large packing, this theory has concluded 
that the spatial distribution of stresses is such that it will 
minimize H within that fabric. 

Now we shall consider the physical meaning of H. By 
its definition, minimum H corresponds to maximum con- 
tact force entropy. Thus, the maximum entropy condi- 
tion has been proven. This proof depends only on the 
stosszahlansatz, which shall be validated by comparing 
its predictions against numerical data. The conclusion 
is that Edwards' hypothesis is not necessary to assert 
maximum contact force entropy. 



IV. A NEW DERIVATION OF THE DENSITY 
OF STATES AND P{f) WITHOUT EDWARDS' 
HYPOTHESIS 

The foregoing granular 7J-theorem tells us that for in- 
finitely large packings maximum entropy (d/da;)-ff = 
persists in every layer, and by extension to very large, 
finite packings it persists in the majority of layers away 
from the boundaries but with greater fluctuations occur- 
ring in smaller packings. The sufficient and necessary 
condition for maximum entropy is given by Eq. 1281 This 
can be written in the form, 



n 



p{ga) 



n 



pig'J 



^^i9») t^i^Hg'o) 



(32) 



which implies that either side of the equation may be 
written as equal to a constant. Taking the logarithm. 



51 log 



pjga) 

T±(5a) 



c 



(33) 



Since this is valid for every possible "rn-nary" transition 
of grains during a translation distance of Aa;, its most 
general solution is when C is written as a linear combi- 
nation of all quantities that are conserved in the trans- 
lation. Since we did not specify the orientation or the 
direction of the translation, this result is valid for all. As 
discussed in Sec.IIA, the total force in every orientation 
must be conserved. Fabric is also conserved by definition 
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of the problem. Therefore, 

log|^ ^4) (34) 

for every € S. The first two terms are for force con- 
servation and the last is for fabric conservation. (We 
can write the first two terms using stress tensor notation 
if we wish, and hence explicitly include shear stresses.) 
Rearranging to obtain p, 

p{g) = G{e,, 0,)^{g)T{g)e-P^^^-^y^y (35) 

where we recall that 5* is the function that enforces 
the bounds on iS, evaluating either to unity or zero if 
the cohesionless grain is stable or unstable, respectively. 
G{9i, . . . ,64) is the fabric partition factor. Eq. [551 is 
identical to the equation derived in Q by Gibbs' most 
probable distribution method. That method assumed the 
maximum entropy condition as a subset of the claims in 
Edwards' hypothesis, but here it has been derived. Re- 
membering that T is a functional of P{f) by Eg . [TT] and 
that 

P{f)= I -Dg p{g) 5\f - f{g)) (36) 
Js 

we see that Eqs. [35] and [36] form a recursion so that p may 
be solved numerically when stress and fabric are speci- 
fied. Solving for p{g) provides everything that can be 
known about the single particle density of states, includ- 
ing P{f), so it makes testable predictions. 

V. EMPIRICAL VALIDATION OF THE 
GRANULAR STOSSZAHLANSATZ 

As with Boltzmann's theory, the translation equation 
and the granular H theorem depend upon the validity of 
the stosszahlansatz. It can be tested only by comparison 
with experimental or simulation data. First we must ob- 
tain predictions from the theory so we have something to 
compare. The numerical solution of p{g) from the recur- 
sion equation has been obtained (in approximation, for 
the isotropic case) and the details of the method are pro- 
vided in 0. The method uses Monte Carol integration, 
randomly sampling the region of integration with a flat 
measure. This is performed in a recursion, with p{g) ob- 
tained from P{f) and vice-versa until convergence. This 
recursion begins from an arbitrary initial condition, for 
example P{f) = 5{f — 1). Several different initial condi- 
tions were checked and all converged to the same result. 
This method was used to produce a representation of 
p{g) consisting of 12 billion grain configurations, which 
was sufficient for very smooth statistics. On an aver- 
age desktop workstation the algorithm converged for a 
sample of one million grains in about one minute, with 
a larger number of grains taking proportionately longer. 
Several different algorithms and different approximations 
to the math were used and they resulted in only minor 



variations in the statistical results, demonstrating the ro- 
bustness of the basic form of the solution. Since a pack- 
ing's fabric P4g{0i, . . . , 64) will evolve as the packing is 
sheared, it is desirable to force the numerical solution to 
a p(g) that has some particular fabric that is found in 
a real case. Thus, it is necessary to weight the Monte 
Carlo sampling to include some classes of contact angle 
configurations more often than others. The weighting 
may be determined iteratively by adjusting the weight 
factor used in the algorithm until the desired fabric is 
obtained. This weight factor, multiplied by zero in the 
regions of steric exclusion, is in fact the fabric partition 
factor G{6i, . . . ,64) that appears in Eq. [35] Similarly, 
the overall stress in the solution p{g) may be driven to 
any desired state by weighting the Monte Carlo sampling 
with the Boltzman factor shown in Eq. 1351 

This numerical solution has been compared to numer- 
ical data from discrete element modeling with idealiza- 
tions approaching those of the theory, and a subset of 
the results have appeared in a Letter [7[ with an archival- 
length paper to follow. The simulation used 17,000 grains 
that were 2D, round, frictionless, and cohesionless. The 
fabric used in the theory implies that they are monodis- 
perse (in that the steric exclusion angle was assumed to 
be precisely 7r/3 radians). However, the simulation used 
a small polydispersity of 1.5 to avoid crystallization of the 
grains. To approach the isotropic idealization the grains 
were deposited into a rectangular, rigid-walled container 
randomly and without gravity, and then expanded in di- 
ameter while allowed to push each other around until 
they jammed. The grains within four grain diameters 
of the walls were discarded from the statistical analy- 
sis to avoid boundary effects, since the theory describes 
the bulk of infinitely large packings (where boundaries 
do not exist). Residual kinetic energy was viscously dis- 
sipated until the forces had negligible dynamical ffuctu- 
ation, representing the idealization of a perfectly static 
packing. The grains had a linear spring contact law, but 
the packing was kept as close as could be achieved near 
the limit of jamming so that minimal compression of the 
contacts occurred. This approximates the perfect grain 
rigidity and the isostaticity of the theory wherein the 
exact form of the contact law becomes irrelevant. The 
data show that with these idealizations the predictions 
of the theory are validated. An example of the corre- 
spondence between theory and empirical results is shown 
in Fig. [6] which shows P{f) obtained from the theory 
and from the DEM data for the Z — 4 grains. Discus- 
sion of the results of Z ^ 4 grains is found in ^ and 
will be expanded in the future archival-length paper. It 
should be noted that there are no free parameters in the 
theory that could have been adjusted to obtain the cor- 
respondence; the good agreement occurs automatically. 
Thus, the stosszahlansatz has been validated along with 
the maximum entropy condition that it predicts for this 
isotropic case. We may conclude that the reductionist ex- 
planation for granular contact force statistics is that they 
exist at maximum entropy, and that is because correla- 
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FIG. 6: Semi logarithmic P (/) obtained in numerical solu- 
tion of the translation equation (theory) and from the Z = 4 
population of a DEM simulation. The statistical fluctuations 
in the tail of the DEM data (amplified on the log axis), do 
not exist in the theory because of the thermodynamic limit. 

tions cannot arise through the closure of loops of grains. 
This statement should be further tested in non-isotropic 
and less ideahzed cases. 



VI. SUMMARY AND CONCLUSIONS 

By following Boltzniann's method but with symmetric 
information flow in the collision network, this paper has 
shown that the layer-by-layer relationships in granular 
contact forces, like the time-evolution of momenta in a di- 
lute gas, must cause their statistical distribution to relax 
to the form that represents maximum entropy. This pa- 
per also derived a version of Boltzmann's stosszahlansatz 
and a version of Shannon's entropy (negative H) relevant 
to granular contact forces, and found a way to derive p{g) 
without recourse to Edwards' hypothesis. It has been 
demonstrated previously with more results to appear in 
a future archival paper that the predictions of p{g) are 
in outstanding agreement with numerical simulations. 

Boltzmann's H theorem for the dilute gas tells us 
that every tiny segment of the Poincare cycle is dom- 
inated by the Maxwellian distribution, so that only a 
tiny amount of dynamics is needed to depart from any 
non-Maxwellian state and settle into maximum entropy. 
Both the Poincare ergodicity and the Boltzmannian re- 
laxation are thus attained by traveling along the tra- 
jectory in phase space. Static granular packings, on 
the other hand, do not travel along any trajectory in 
phase space. The granular translation theorem devel- 
oped above depends upon contact forces maintaining 
static spatial relationships, not traveling and interact- 
ing through time. Whereas Boltzmann's proof obviated 
the need to travel the entire Poincare cycle and showed 
that the system need travel only a tiny segment of it 
to justify the assumptions of statistical mechanics, the 
granular H theorem shows that a granular packing need 
not travel through F space at all, because it has special 



relationships between its own layer-wise subspaces built 
into its single locus in F space. The author proposes 
that this special feature of granular packings deserves 
a name, and "self-ergodic" seems appropriate. That is, 
the system must exist in a relaxed state by nature of its 
internal relationships (self-enforced and affecting itself), 
and this produces the same statistical characteristics that 
an ergodic theorem seeks to establish for kinetic systems 
exploring all conserved-energy states with equal proba- 
bility. Thus the terminology: self -I- ergodic. 

There are of course many significant classes of granular 
packings that will not be relaxed to maximum entropy, 
but these have not been discussed here. These include 
the boundary regions near the container walls of a pack- 
ing, as well as granular packings that were prepared to 
have abrupt changes in fabric somewhere within their 
bulk so that the density of states cannot be in its most 
relaxed state within the narrow band of grains on both 
sides of the interface. The difhculty in applying this the- 
ory to boundary regions is that the fabric does not remain 
constant near the boundary, whereas the theory assumes 
that the fabric is constant. Another case not described 
by the theory is the thin layer of grains at the top surface 
of a packing in gravity, where the self-weight introduces 
non-negligible stress gradients. However, solution of the 
more symmetric case considered here opens the door to 
solving more complex problems. 

This theory has been developed only for the case of 
rigid, round grains without friction. This avoids the com- 
plication of rotational degrees of freedom for the grains. 
The author does not know how to define the grain states 
to include torques and rotational degrees of freedom so 
that they would have sufficient symmetry for the anal- 
ysis. The assumptions of the theory also exclude cases 
with highly compressed grains, which are far from iso- 
staticity [15]. 

The theory has been developed assuming Z ^ A for 
every grain, although frictionless 2D packings with low 
polydispersity are known to have grains with Z — 3, 
Z = 5 and occasionally Z = 6. It seems possible to 
generalize the theory to allow for different values of Z 
as long as the average is still 4 for isostaticity. In this 
paper, the grain state g was defined to have just the two 
variables with units of force, Wx and Wy, plus the four 
contact angles. We may alternatively use the two vari- 
ables t = Wx + Wy and s — {wx — Wy)/t. So, for grains 
with z = 3 we would define gs — {t, 9i,92, ^3}, and this 
would be sufficient to define the grain's entire state al- 
lowing us to solve all the contact force values on the grain 
by linear algebra and trigonometry. For Z = 4 we would 
use (74 = {i, s • i, 01, . . . , ^4}, which is equivalent to the 
treatment given in this paper. Then, for Z — b we need 
one additional state variable in (75 to allow us to solve all 
five contact forces on the grain. This additional variable 
could simply be the value of one of the contact forces, 
but to maintain symmetry among the contacts it would 
presumably be better to define something analogous to 
a quadrupole term and thus continue the sequence of t 
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and st, which represent the monopole and dipole terms, 
respectively. Obtaining a numerical solution to the re- 
sulting theory would be a simple extension of the existing 
Z = 4 algorithm. 

It is unknown whether this theory will work for an or- 
dered (crystalline) granular packing. Experimental and 
numerical studies have shown that even the microscopic 
variations in packing geometry of crystalline packings are 
sufhcient to break the symmetry and relax P{f) to a typ- 
ical form . This theory may take advantage of this by 
keeping translation distances Ax small enough that the 
microscopic packing variations are orders-of-magnitude 
larger. Thus even with a crystalline packing the assump- 
tions behind Eqs . fTOl and fT2l should be valid and the the- 
ory may work. Nevertheless, it is possible that corre- 
lations will arise in the regular pattern of closed loops, 
violating the generalized stosszahlansatz. It would be in- 
teresting to study how much disorder is required for the 
theory to work. 

One application of this theory in particular was hy- 
pothesized by the author several years ago [l3| . The idea 
is that, since granular contact forces are in a state of equi- 
librium at maximum entropy subject to the layer-by-layer 
conservation of forces in each direction, then it is possible 
to define a rank-2 tensor "temperature" for the contact 
forces and even a coordination number "chemical poten- 
tial" to explain the partition of stress and fabric fluc- 
tuations throughout a granular packing. This approach 
may lead to a full theory of rheology. These observations 
were apparent when numerical solutions to the theory 
first proved robust and convergent, even before this for- 
mal proof of maximum entropy was accomplished. That 
is because the numerical convergence discussed in the 
earlier publications was in fact an empirical demonstra- 
tion of the same concept. Once the "self-ergodic" nature 
of granular materials is identified, then the temperature 
and entropy concepts fall out rather straightforwardly. 
A future publication will be forthcoming to explain these 
concepts along with a series of discrete element modeling 
simulations that have been performed to test them and 
to draw further conclusions. 



APPENDIX A: A PRIORI ARGUMENTS FOR 
THE GRANULAR STOSSZAHLANSATZ 

As illustrated in Fig. [71 correlations between neighbor- 
ing contacts on the same grain arise either through the 
grain itself or through the loops in the packing. A typ- 
ical loop is four or more grains, so going the long way 
around a loop induces a relatively weak four-point corre- 
lation, but going the short way between the two contacts 
(staying intra-grain) induces a much stronger two-point 
correlation. 

Furthermore, we may note that small loops composed 
of N grains, 3 < < 8, are formed through adjacent 
pairs of contacts that are closer to the non-correlating 
7r/2 angle than the correlating tt radians of separation 




FIG. 7: Correlations between neighboring contact forces may 
arise two different ways. (Left) They arise through intra-grain 
stability requirements, resulting in a strong, two-point corre- 
lation between the forces at B and A. (Right) They arise 
through a series of intra-grain stability requirements work- 
ing grain- by- grain around the loops, resulting in weak higher- 
order correlations. The example here shows a four-point cor- 
relation through the loop, _B to C to D to A, which will be 
much weaker than the direct two-point correlation from B to 
A. 




FIG. 8: Contacts that are highly correlated are close to tt 
radians apart on the same grain [^. Hence, a closed loop 
composed of highly correlated pairs of contacts can turn only 
very slowly and must pass through a very large number of 
grains. 

(see 0), and hence there should be minimal correlation 
in each two-point leg of a force loop. The composite cor- 
relation going all the way around the loop must there- 
fore be exceedingly weak. On the other hand, force loops 
that pass through a larger number of grains and hence 
closer to tt radians of separation for each two-point leg 
of the loop will require a vastly larger N to slowly turn 
through the full 27r radians to close the loop as illustrated 
in Fig. [SI Because N is so large {N > 8), we would there- 
fore expect these A^-point correlations to be very weak, 
as well. 

The total correlation between a pair of contacts on a 
grain must be the sum of information from all the loops 
in the packing that contain the grain in question. Due 
to the disorder of the packing and the large number of 
loops that contain the same grain, some correlating and 
some anti-correlating, it is expected that the contribu- 
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tions from increasingly larger loops of grains will be in- 
creasingly decoherent and largely cancel one another. 

From these arguments there is good reason to as- 
sume a priori that the intra-grain contribution to the 
correlations is the dominant one and that closed-loop 
contributions may be discarded. This will be the 
granular stosszahlansatz. The ultimate proof of this 
stosszahlansatz is its ability to make valid predictions. 
Comparison with empirical data 0, 0] have validated its 
predictions. 



APPENDIX B: TWO FORMS OF THE 
COLLISION AUXILIARY IN THE DENSITY OF 
STATES 

The purpose of the "collision auxiliary" term in Eq.[¥]is 
to restrict the allowable set of vectors "colliding" (meet- 
ing together) at the nodes (grains or particle collisions). 
Already the vector sums of the forces or momenta are 
automatically conserved (a la Newton's Second Law) by 
the proper selection of phase space coordinates. However, 
in both granular mechanics and in the kinetic theory of 
dilute gases there are additional restrictions imposed by 
physics upon these sets of vectors. 

For the case of the dilute gas, in addition to conser- 
vation of momentum there is the requirement for the 
conservation of energy, and even that is not sufficient 
to define the allowable sets of vectors appearing at each 
collision. A pair of colliding momentum vectors may 
meet at any arbitrary magnitudes and angle relative to 
one another. The magnitudes and directions of the two 
outgoing momentum vectors are then determined by the 
precise form of the interaction potential between the par- 
ticles. Thus, the collision auxiliary term in Eq. 3] must 
be written to specify the set of allowable outgoing vec- 
tors based upon the incoming vectors and the interaction 
potential. 

For the case of granular packings, on the other hand, 
the sets of allowable force vectors appearing on each grain 
are not so precisely constrained. In addition to the con- 
servation of force from one hemisphere to the next (New- 
ton's Second Law) there is only the requirements of steric 
exclusion. There is no analog of interaction potential or 
conservation of energy to determine two of the vectors 
as a function of the other two. Precisely because of this 
additional freedom, granular packings exhibit memory in 
their fabric, the statistical preponderance of their con- 
tact angles resulting from past disturbances of the pack- 
ing. Therefore, an ensemble of granular packings would 
need to have the current state of its fabric specified in or- 
der to be a completely defined ensemble. The "collision 
auxiliary" is therefore used for this purpose. 

Although the collision auxiliary plays two very differ- 
ent roles in granular versus kinetic ensembles, they have 
been lumped together here under the name "collision 
auxiliary" because they both play the role of defining the 
sets of vectors appearing at the nodes insofar as required 



by the physics. The important point to note is that it 
does not matter which form we use in this paper, either 
the kinetic or the granular form, because this term factors 
out from the sum in Eq.[S]and does not affect the results. 
Thus, the results we obtain are valid for either granular 
or kinetic ensembles and we can compare the mathemat- 
ical form derived for the stosszahlansatz in one case with 
the mathematical form derived in the other. 



APPENDIX C: DERIVATION OF THE 
STOSSZAHLANSATZ FOR TWO CASES OF 
INFORMATION FLOW 

Beginning from Eq. [71 we expand the product over Ci, 

Nj\ Nj 2 



hm nEnn 

Nj 4 

X n Il'^c. s^r-'^fP') (CI) 

[3j—aj (5—3 



where the Kronecker delta function is nonzero only for 
contact pairs contained within list C;. We have ordered 
the contact pairs in Ci so that the grain in the direction of 
lower j is a and contacts in the increasing j hemisphere 
are numbered 1 and 2. Note that the product in /? is only 
from a to N, the upper triangle in the (a, /3) plane, since 
Newton's third law is enforced only once for every pair, 
and this produces the correct number of delta functions 
so that the units of p are correct. In writing it this way, 
we have preferentially given it an asymmetry in the j 
direction. 

Now we must make the critical assumption about 
causality in this network. First we attempt to recover 
Boltzmann's stosszahlansatz for the dilute gas. We do 
so by treating the information asymmetrically in the j 
direction, so that the probability of finding a certain set 
of grains in j depends only on the set of grains in j — 1 . 
Therefore, for this case we may start from Eg. lCll Com- 
muting the sum in i with the first two products, 

$ = lim IT TT TT y 

Nj 4 

X n I[sc,s^ir'' + f') 

f3j—aj 5—3 

Nj 2 Nji. 

= hm A^nnriE 

Nj 4 

X E E ^'(r^ + O (C2) 

f3j—aj (5—3 

This was renormalized by M since each of the terms in 
the product is now a sum of many delta functions instead 
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of just one delta function. This also introduced a small 
error inside the limit in that the product of the sums 
expands to a sum of the product of all possible combina- 
tions of the delta functions with replacement, whereas it 
should have been without replacement. Taking the limit, 
the cross-terms with non-representative distributions of 
single grain states pg become a set of zero measure and 
so the error vanishes. Next, we expand the delta function 

N, 2 Ni\ 



$ = hm X TT TT TT 



j aj = l7=lij=l 



Ni 



X lim - — lim — — 

^ ^ Af^o Af Ae-,0 A9 

f3j=aj (5=3 

X {e - fkin)] - e - h+iin)] } 
X {e [|6"^'^-6',(r'')| -tt] 

-e[|0''^-%i(r^)|-7r]} (C3) 

where the subscripts k and / define the k^^ interval of 
the force axis as the interval containing /"''' and the l^^ 
interval of the angle axis as the interval containing 6°''^. 
Taking the summations. 



TV 



M TT lim lim lim -i- -r;: TT TT ' 

-"■-"-{Jva^loo} A/^oAe^o AF A9 

j a—l 7=1 



,0") 

'kl 



(C4) 

where n^j^i is the number of contacts in region j that 
fall into the /c"^ and Z'^ bins. Taking the limits converts 



(C5) 



This says that the probability of finding a member of the 
ensemble with a particular set of force pairs intersecting 
on a common grain in the j*"^ layer is Ho'^i Pjif^^ : 
where the distribution of intersecting forces in the j^^ 
layer of the ensemble is. 



We have derived Boltzmann's stosszahlansatz. We ob- 
tained it from a statement about the topology of the 
network — no closure of force (or momentum) loops — 
rather than simply writing it as a statement borrowed 
from probability theory. 

Now we can put into the analysis the symmetry that 
is correct for a granular packing, so that causality op- 
erates equally in all directions. We do so by modifying 
Eq. ICll by including all contacts on each grain and by 
including the entire (a, [3) plane, so that every contact is 
considered twice: once looking in the +j direction and 
again looking in the — j direction. (We can also relax the 
restriction that contacts 1 and 2 must be on a particu- 
lar hemisphere of the grain with contacts 3 and 4 on the 
opposite hemisphere. This more symmetric treatment 
allows grains to have three contacts on one side of the 
grain and only one contact on the opposite side.) But 
now that we are multiplying over the entire (a, /?) plane 
there are twice as many delta functions and so in effect 
$ has been squared. 



Nj>. Nj 4 

$2 ^ n E n n 



Nj 4 
I3j = l S=l 



(C7) 



(C6) 



The analysis is now symmetric in j and proceeds identi- 
cally as before, with the end result, 

$2 = X n n PAr')PAr^)PAr')W) (cs) 

j a=l 

SO that by taking the square root of both sides we con- 
clude 

F,{w^,wy,6^, ... ,04) = [P-Afi) PAh) Pjifs) PAh)]"^ 

(C9) 

This is the granular stosszahlansatz. Both this version 
and Eq. IC6I say that the topology of the network does 
not pre-correlate the intersecting force vectors. 
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